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I . PROBLEM  STATEMENT 

1.1  Backg round  - The  unsteady  exhaust  plume  and  blast  field  genera- 
ted by  a rocket  launch  has  been  of  concern  in  the  past  principally  to  those 
responsible  for  the  structural  integrity  of  the  launch  pad  and/or  adjacent 
structures.  The  transient  thermal  loads  and  forces  may  be  conservatively 
estimated  without  recourse  to  detailed  analysis  of  the  plume  flowfield.  How- 
ever, the  steady  exhaust  plume  flowfield  has  been  studied  in  careful  detail 

at  both  high  and  low  altitude,  since  it  may  provide  a trail  of  optically  and/or 
electromagnetical ly  observable  gases  suitable  for  remote  detection  and  tracking 
of  a launch  vehicle,  and  possibly  diagnostic  information  as  well.  The  objective 
of  early  detection  has  pressed  interest  to  the  lowest  practical  altitudes,  but 
the  unsteady  exhaust  plume  and  blast  field  at  ignition  has  not  been  addressed 
in  this  context  since  any  signal  (optical,  electromagnetic  or  acoustic)  is 
largely  attenuated  by  the  atmosphere  near  sea  level.  However,  the  plume  tech- 
nology base  is  nevertheless  incomplete  without  an  adequate  description  of  the 
plume  flowfield  at  ignition  and  lift-off. 

1.2  General  Description  of  the  Flowfield  - The  unsteady  flowfield  at 
ignition  of  a rocket,  as  depicted  in  Figure  (1)  bears  certain  resemblances  to 
the  muzzle  blast  field  of  a gun  or  mortar,  shown  in  Figure  (2).  Both  share 
structural  features  common  to  the  steady  underexpanded  rocket  exhaust  plume; 
principally,  the  now-familiar  underexpanded  supersonic  core,  barrel  shocks  and 
Mach  disc.  The  steady  exhaust  plume  is  most  frequently  depicted  in  a vehicle- 
fixed  frame  of  reference,  with  the  vehicle  far  away  from  the  ground  surface.  In 
this  situation,  the  exhaust  plume  near  the  vehicle  is,  on  the  average,  steady 
and  predominately  supersonic.  Turbulent  mixing  eventually  decelerates  the  plume 
to  entirely  subsonic  conditions  further  from  the  vehicle,  and  asymptotically  far 
from  the  vehicle  to  a randomly  fluctuating  field  to  be  dispersed  by  prevailing 
winds,  buoyant  forces  or  simply  dissipated  by  viscous  stresses.  However,  the 
picture  is  considerably  altered  at  launch,  due  to  the  interaction  of  the  plume 
with  the  ground  plane.  In  simplest  terms,  an  unsteady  plume  spreads  over  the 
ground  surface,  displacing  the  air  around  it.  The  unsteady  plume  collapses  as 
the  vehicle  leaves  the  vicinity  of  th<»  ground  plane,  and  the  remaining  propel lant 
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FIGURE  1.  STRUCTURE  OF  IGNITION  PHASE  BLAST  FIELD 
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FIGURE  2.  STRUCTURE  OF  BLAST  FIELD  OF  A SMALL 
CALIBER  (5.56mm)  RIFLE 
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gas  cloud  is  eventually  dispersed  by  prevailing  winds  and  buoyant  forces. 

During  the  interim  period,  however,  the  unsteady  plume  will  produce  a blast 
field  and  flash  analogous  to  (although  perhaps  less  intensely  than)  the  well- 
known  muzzle  blast  and  flash  of  a gun.  The  strength  of  the  blast  field,  as 
measured  in  terms  of  overpressure,  will  depend  on  factors  like  the  ignition 
time,  vehicle  acceleration,  etc.  The  flash  intensity  will  depend  on  these 
factors  plus  afterburning  of  unburnt  propellant  and/or  combustion  products, 
opacity  of  the  surrounding  layers,  etc. 

1.3  Present  Approach  - Obviously,  the  subject  problem  is  very  complex 
to  describe  in  complete  detail.  However,  it  can  be  reduced  to  an  assembly  of 
component  problems  for  which  approximate  solutions  are  tractable,  which  yield 
physical  insights  to  the  overall  problem,  and  from  which  limited  scaling  laws 
suitable  for  certain  applications  can  be  derived.  The  present  report  summarizes 
the  results  of  study  of  the  various  component  problems  and  attempts  to  bring 
them  together  into  a un i f ied,  although  approximate,  description  of  the  main  fea- 
tures of  the  unsteady  plume  and  blast  field  generated  by  a missile  launch. 

The  present  approach  to  the  subject  problem  derives  in  a large  part  from  previous 
work  or,  muzzle  blast  field  analyses  for  guns^>2.  These,  in  turn,  initially  view- 
ed the  muzzle  blast  as  a spherical  explosion,  with  a constant  rate  of  energy  aJ- 
dition^.  This  approach  successfully  described  the  trajectory  of  the  leading 
shock  wave  at  early  times,  and  thus  lent  credibility  to  further  exploitation  of 
the  assumption  of  spherical  symmetry.  The  asymptotic  shock  motion  was  subse- 
quently^  described  by  Whitham's  theory**.  This  approach  also  clearly  identified 
the  rate  and  duration  of  energy  addition  as  primary  parameters  controlling  the 
blast  field.  Accordingly,  the  rate  and  duration  of  energy  deposition  at  the 
ground  surface  from  an  accelerating  rocket  Is  derived  In  Section  III.  The  asymp- 
totic inviscid  structure  of  the  plume  and  blast  field  is  reviewed  In  Section  IV. 
Description  of  the  actual  thermochemical  properties  in  the  unsteady  inviscid 
shock  layer  requires  numerical  solution  of  the  governing  equations;  the  method- 
ology and  sample  solutions  are  presented  in  Section  V.  A simple  model  of  the  un- 
steady turbulent  mixing  layer  that  develops  along  the  Interface  between  the  ex- 
haust gases  and  the  air  has  been  developed  and  is  discussed  in  Section  VI. 
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Finally,  in  Section  Vll  an  overview  of  the  problem  is  presented  in  terms  of  the 
model  derived  herein,  and  conclusions  are  discussed. 
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1 II.  DESCRIPTION  OF  THE  PRESENT  MODEL 


2.1  Time  Scales  - A rocket  motor  begins  to  transfer  thermal  and 
kinetic  energy  to  atmosphere  upon  emergence  of  the  leading  pressure  wave  from 
nozzle.  After  a series  of  reflections  of  pressure  waves  within  the  combustion 
chamber,  an  approximately  constant  chamber  pressure  and  burning  rate  is  achieved, 
if  the  combustion  process  is  stable,  producing  correspondingly  constant  thrust 
and  energy  flux  rate  in  the  exhaust.  The  ignition  transient  period  (for  a 
solid  propellant)*  is  of  the  order®: 

tj/C|,/L*  - 10  to  20  (1) 

For  example,  using  L*  = 10  feet  and  c^  “ 3000  ft/sec  as  representative  values, 
give  ignition  times  of  about  30  to  60  milliseconds. 
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As  the  vehicle  lifts  off,  the  rate  of  energy  deposition  in  the  atmosphere  at 
the  launch  site  varies  due  to  a combination  of  relative  motion  between  the 
vehicle  and  the  ground  surface,  and  turbulent  mixing  of  the  exhaust  gases  and 
the  surrounding  air.  For  the  present  purposes,  lift-off  is  considered  complete 
when  the  rate  of  energy  deposition  at  the  launch  site  drops  to  zero.  The  lift- 
off period,  tj^,  is  thus  defined  as  the  time  from  ignition  to  completion  of  lift- 
off. It  will  be  shown  subsequently  that  this  period  ranges  from  about  5 seconds 
for  a high  acceleration  (100  g)  launch  to  about  12  seconds  for  a low  accelera- 
tion (1.2  g)  launch.  In  the  high  acceleration  case,  t^^  'v  a S thus  even  a 
five-fold  increase  in  acceleration  (500  g)  will  still  result  in  t^^  > lOtj. 
Therefore,  for  a very  wide  range  of  conditions  of  interest,  ignition  is  virtually 
instanteous  on  the  time  scale  of  lift-off,  and  the  two  proecesses  may  be  con- 
sidered independently. 

2.2  Flow  Pattern  - In  view  of  the  different  time  scales  for  ignition 
and  lift-off,  the  interaction  of  rocket  plume  and  the  ground  plane  can  be 
regarded  as  the  impingement  of  an  established,  steady  plume  on  a moving,  planar 


Liquid  propellant  ignition  transient  periods  have  not  been  examined,  but  should 
be  of  a comparable  order,  i.e.,  the  time  for  about  10  to  20  wave  reflections  in 
the  combustion  chamber. 
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^ surface,  as  indicated  schematically  in  Figure  (3).  The  flow  is  turned  from 

i 

vertical  to  horizontal  through  a shock  (Mach  disc)  and  subsonic  stagnation 
point  type  fiow  pattern.  The  flow  then  expands  radially  outward  along  the 
ground  surface.  The  flow  can  either  continue  to  decelerate,  thus  remaining 
subsonic,  or  can  re-expand  becoming  supersonic  very  quickly;  the  former  case 
leads  to  a continuously  increasing  pressure  which  in  general  cannot  match  the 
prevailing  atmospheric  pressure,  whereas  the  latter  gives  rise  to  a cylindri- 
cally  symmetric  supersonic  plume  pattern  (as  indicated  in  Figure  1)  in  which 
the  flow  shocks  up  to  the  prevailing  back  pressure,  through  a Mach  disc,  as  in 
the  case  of  a muzzle  blast  field^.  Thus  the  supersonic  plume  pattern  is  as- 
sumed to  prevail,  in  general.  This  pattern  has  also  been  observed  experiment- 
ally in  this  case  steady,  underexpanded,  sonic  jet  impinging  on  a wall®. 
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Since  the  ignition  time  is  short  relative  to  the  lift-off  time,  the  flowfield 
can  be  considered  to  begin  impulsively.  Accordingly,  a leading  pressure  wave 
Is  driven  ahead  of  the  expanding  propellant  gases  into  the  ambient  atmosphere. 
Under  suitable  conditions,  this  wave  can  be  represented  by  a semi-spherical 
blast  wave,  which  asymptotically  deteriorates  to  a sound  wave,  having  it  origin 
at  the  point  of  impingement  of  the  exhaust  gases  on  the  ground  surface.  As  in 
the  muzzle  blast  case  , the  effective  origin  may  actually  be  a ring  within  the 
cyl indr ical ly  symmetric  supersonic  plume  which  forms  along  the  surface,  lead- 
ing to  a distortion  of  the  assumed  semi-spherical  shock  shape,  but,  like  the 
muzzle  blast,  the  distortion  is  probably  most  pronounced  at  very  early  times 
and  asymptotically  negligible. 

As  in  a muzzle  blast,  the  flowfield  is  quickly  engulfed  by  turbulence,  and  the 
inviscid  structure  shown  in  Figure  (3)  at  time  t^  is  covered  by  a turbulent 
cloud  as  indicated  at  time  t2>  In  the  present  problem,  the  exhaust  plume  shear 
layers  present  one  evident  source  of  turbulence  due  to  velocity  jumps  in  the  in- 
viscid flow.  However,  as  in  the  muzzle  blast,  the  extreme  temperature  jump 
across  the  radially  expanding  front  between  the  propellant  gases  and  air  (termed 
the  contact  surface  in  Figure  1)  presents  another  source  of  turbulence.  Although 
the  effects  of  turbulent  mixing  are  important,  the  inviscid  structure  of  the 
flowfield  must  be  determined  first.  The  turbulent  mixing  processes  can  then  be 
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FIGURE  3.  SCHEMATIC  OF  INTERACTION  OF  STEADY  EXHAUST  PLUME  WITH  THE 
GROUND  SURFACE 
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superimposed  on  the  invlscld  solution,  In  the  usual  vein  of  boundary  layer 
theory,  as  a first  approximation. 
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^ III.  ENERGY  DEPOSITION 

3.1  General  Remarks  - In  an  earth-fixed  frame  of  reference,  the 
rocket  motor  can  be  regarded  simply  as  a source  of  high  temperature,  high 
pressure  gases  released  into  the  atmosphere  in  an  confined  stream  which  im- 
pinges on  the  ground.  The  point  of  impingement  can,  in  turn,  be  viewed  as  an 
effective  source  of  energetic  gases  which  spread  over  the  ground  as  an  un- 
steady plume  and  blast  field.  This  source  is  fed  by  the  exhaust  stream  of  the 
rocket  motor.  Since,  as  noted,  the  ignition  transient  time  is  short  relative 
to  the  lift-off  time,  the  process  can  be  considered  to  begin  impulsively. 


In  this  source-flow  approximation,  the  similitude  solutions  of  blast  wave  theory 
apply,  and  the  boundary  conditions  can  be  stated  in  terms  of  one-dimensional 
flow  representations  of  the  rocket  exhaust  gases.  The  boundary  conditions  con- 
trolling the  supersonic  plume  can  be  approximated  from  the  same  one-dimensional 
model.  However,  the  blast  wave  solutions  can  be  stated  in  terms  of  the  rate  of 
energy  deposlton,  whereas  the  plume  additionally  requires  the  knowledge  of  the 
mass  flow  rate,  since  the  plume  Is  two-dimensional  while  the  blast  solution  is 
one-dimensional . 

3.2  Inviscid  Seal ing  - A brief  review  of  the  essential  elements  of 
the  theory  of  variable  energy  blast  waves3  is  informative  with  respect  to  the 
subject  problem.  In  particular,  the  energy  released  at  the  source  of  the  blast 
wave  must  conform  to  the  power- law; 

E(t)  - Wt®  (2) 
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n - (6+2)/(j+3)  (5) 

Jq  (n.j.y)  ” constant,  and  B » constant  (6) 

The  Index  j Is  2 for  spherical  symmetry  (as  In  the  present  case),  1 for  cy- 
lindrical symmetry  and  0 for  planar  symmetry.  Thus,  with  j » 2 and  3*0  the 
classical  point  explosion  solution,  n * 2/5,  Is  recovered,  whereas  with  6 * 1 
the  solution  for  a constant  rate  of  energy  addition,  n = 3/5,  Is  obtained. 

Few  actual  processes  of  interest  conform  to  a power  law  variation,  although  a 
combination  of  different  power  laws  for  different  periods  can  provide  a reason- 
able approximation  to  the  total  energy  release  rate.  In  the  present  problem 
the  combination  of  3 = 1 for  0 ^ t ^ t^  8 < 1 for  t ^ t2  and  3 = 0 for 

t > t2  will  prove  to  be  appropriate.  However,  the  validity  of  blast  wave 
theory,  or  of  scaling  based  on  blast  theory  variables,  may  be  questioned  if 
8 ^ constant.  Oabora^  suggested  use  of  a length  scale  based  on  the  total 
amount  of  energy  released  as  t In  the  present  case  this  length  scale  Is 

R^  - (2E(t2)/P„)’^^ 


where  dE/dt  * 0 for  t > t2,  and  the  factor  of  2 has  been  Inserted  to  account 
for  symmetry  with  respect  to  the  ground  plane.  Finite  difference  solutions 
for  variable  3 and  widely  differing  values  of  R^  were  examined  In  Reference  (2); 
the  trajectories  of  the  blast  wave,  and  the  contact  surface  were  found  to  corre- 
late according  to: 

r /R  - f,  (t  c /R  ) (8) 

so  I CD  O 


r /R 
c o 


\ (« 


(9) 


for  t < 12*  Examples  are  shown  In  Figures  (4)  and  (5),  where  the  data  points 
are  taken  from  the  numerical  solutions  of  Reference  (2)  for  various  size  and 
type  high  muzzle  velocity  weapons.  The  Mach  disc  pos  it  ion  was  al so  found  to  be 
correlated  by: 
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r /R  - f,  (t  c /R  ) (10) 

mo  3 oo  o 

for  t <<  12*  The  behavior  of  the  Mach  disc  for  t ->  t2  will  be  discussed  later. 

For  values  of  t >>  t2  the  shock  trajectory  decays  to  the  weak  shock  solution  of 
Whi  tham**  : 
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c t + 


log^  (r^/R^) 


c t 


CO 


(11) 


while  the  contact  surface  tends  toward  a constant  of  order  R . However,  the 

o 

Mach  disc  collapses  back  to  origin  as  t ■+  t2*  It  is  pointed  out  parenthethi- 
cally,  that  as  this  limit  is  approached,  the  pressure  becomes  atmospheric 
everywhere  except  immediately  adjacent  to  the  shock,  and  the  (frozen)  tempera- 
ture of  the  propellant  gases  approaches 

T = T exp  (AS/C  ) (12) 

00  V 

where  AS  refers  to  the  entropy  increment  in  the  propellant  gases  produced  by 
the  combination  of  combustion  and  shock  heating.  Thus  a substantial  tempera- 
ture increment  across  the  contact  surface  will  persist  (in  the  absence  of  mix- 
ing ) as  t ®. 

Therefore,  determination  of  the  rate  of  energy  deposition  is  essential  to  speci- 
fication of  the  boundary  conditions  on  the  subject  problem  and  of  the  length 
scale  for  the  inviscid  structure  of  the  flowfield. 


3.3  Mass  Flow  Rate  - The  energy  flow  rate  in  the  rocket  exhaust  is 
given  by  the  product  of  the  total  specific  energy  and  the  mass  flow  rate.  For 
this  reason,  and  another  to  be  discussed,  it  is  of  Interest  to  consider  first 
the  rate  of  mass  flow  impinging  on  the  ground  surface: 


(u-V) rdr 


(13) 
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where  is  the  plume  radius,  u the  gas  velocity  in  vehicle-f ixed  coordinates, 

‘ and  V the  vehicle  velocity,  or  equivalently  the  velocity  at  which  the  ground 

surface  recedes  from  the  vehicle.  In  the  absence  of  mixing,  the  mass  flow  rate 
in  vehicle-fixed  coordinates  is  constant: 


r 

e 


Thus 

I m - <m>  (1-V/u^)  (17) 


The  vehicle  velocity  is  easily  obtained,  assuming  a vertical  launch: 


TR  246 


and 


< m > = Th/g  Isp 


the  inviscid  approximation  to  the  mass  flow  rate  is: 


Th 

g Isp  ^ ' Isp 


[ , + t . ,og  (1  + (^)  )] 


m^g  Isp 


Th  f , t /Th  , t /Th  -1 

g Isp  Isp  mg  ^ Isp  'm  g) 


(21) 


(22) 


(23) 


This  expression  yields  m = Th/gisp  = constant  when  t/Isp  « (Th/m  g-l)  \ and 

Th  t t ^ 

m ->  0 as  [log(l  +{ ) = — ) - ■= — 1 1 ; i .e.  , a very  slow  rate  of  decrease 

^ mg  Isp  Isp 

in  the  mass  flow  ra?e  at  the  ground  surface.  Returning  to  Equation  (17).  it  is 
obvious  that  in  this  inviscid  approximation,  m -►  0 only  as  the  vehicle  velocity 
approaches  the  exhaust  gas  velocity. 


Mixing  can  be  expected  to  have  two  important  but  opposing  effects  on  the  in- 
viscid result  expressed  by  Equation  (17)  or  (22).  First,  the  mass  flow  rate 
in  the  plume  <m>  will  increase  due  to  entrainment  of  additional  mass  from  the 

atmosphere,  and  second , the  (mass-averaged)  exhaust  gas  velocity  u (z)  will  de- 

e 

crease  due  to  momentum  transfer  with  the  atmosphere.  These  effects  have  been 
examined  by  using  actual  mixing  calculations^  for  a full-scale,  single  nozzle, 
rocket  motor  producing  approximately  2 x lO^N  thrust  at  sea  level,  with  an  ef- 
fective exhaust  velocity  of  u = 2800  m/s  and  a nozzle  exit  radius  of  about 
0.9m.  The  conversion  from  vehicle-fixed  (axial  distance, z)  to  ground-fixed 
(time,t)  coordinates  was  made  using  Equation  (l8),  i.e.: 


Vdt  •“  g Isp  t < 


Th 

m g'  Isp 


1 + (:!:^)  I- 
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m^g  Isp 


m g 
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Representative  curves  of  non-dimensional  mass  flow  rate,  i.  e. , m/<m>,  are 

shown  in  Figure  (6),  using  Isp  = 250  s and  Th/m  g = 1.17,  1.25,  1.50  and  2.0C. 

o 

It  is  evident  that  mixing  first  produces  a marked  increase  in  the  mass  flow 
rate  relative  to  the  inviscid  value,  up  to  nearly  a factor  of  10  in  the  lowest 
acceleration  case,  but  then  an  abrupt  cut-off  occurs  as  u^  is  decreased  and 
approaches  V.  As  the  acceleration  rate  increases,  the  entrainment  of  addi- 
tional mass  decreases,  but  the  cut-off  occurs  earlier  and  more  abruptly.  It 
Is  also  evident  that  the  viscous  solution  at  Th/m^g  = 2.0  bears  no  resemblence 
to  the  corresponding  inviscid  solution.  The  range  of  validity  of  the  inviscid 
solution  will  be  explored  further  in  the  following  discussion;  the  point  here 
is  to  emphasize  the  effects  on  the  mass  flow  rate  produced  by  mixing. 

3.1*  Energy  Deposition  Rate  - The  energy  flux  at  the  ground  sur- 
face is  given  by: 
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E(t)  = i 

2it  1 p(u-V) 

Ea,  ( / Cvj  dT  + Ae^)  + J 
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where  the  contributions  of  thermal,  chemical  and  kinetic  energy  have  been  in- 
cluded. The  corresponding  energy  flux  in  the  steady  exhaust  plume,  in  vehicle- 
fixed  coordinates,  is  given  by: 


<E(z)>  - 2it 
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r T 
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/ pu 

^“i  ( 
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= ) + T 
o z 


rdr 


(27) 


Again  assuming  inviscid  flow,  q = u = gisp,  and  e = h <<  iq  , gives: 
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<E(z)>  = <E  > 

o 


P A • 
9 Isp  (:j;Vf  + 2 


(28) 


Thus: 


E(t)  a < E^>  (1  - V/g  Isp) 


(29) 


and 


E(t) 


E dt  < E > t ( 2 + ^ 


2Isp 


(30) 


1 + (I^)  f 

'm^g'  Isp 

(—  ) - 
m^g  Isp 


fs,^ 


In  accord  with  the  notation  adopted  earlier,  dE/dt  -*■  0 for  t > t2.  and  E(t2) 
denotes  the  total  amount  of  energy  deposited  at  the  launch  site.  Then,  from 
Equation  (30),  E(t)  = E(t2)  for 


log  (1  + (^)  ^ ) - ^ - 1 


(31) 


If  t,/Isp  <<  1 and  Th/m  g >>  1 (i.e.,  a hlqh  acceleration  launch),  then: 

* o 


t2  “ (e  - l)  Isp  (Th/m^g)  ' » 1.718  Isp  (Th/m^g) 


-1 


(32) 


and 


E(t2)/<E^>  a Isp  (Th/m  g) 


5.905  (?-)  + 0.718 
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or 


E(t2) 


= g Isp 


P A 

2 - r 


5.905  {—)  + 0.718 
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=0.36  (g  Isp)' 
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Equations  (32)  and  (35)  provide  very  simple  formulas  for  the  time  required  to 
complete  the  deposition  of  energy  at  the  1 aunch  s i te,  i.e.,  thel  ift-off  time,  and 
the  total  amount  of  energy  deposited,  which  may  be  compared  to  the  "exact"  in- 
viscid  solution  given  by  Equation  (30) . Solutions  to  Equation  (30)  for  a 
range  of  values  of  Th/m^g  are  presented  in  Figure  (7).  The  values  of  t2  and 
E(t2)  given  by  Equations  (32)  and  (35)  are  also  shown  in  Figure  (7)  as  data 
points  for  values  of  (Th/m^g-l)  = 100,  50,  20  and  20.  It  can  be  seen  that 
Equations  (32)  and  (35)  predict  t_  and  E(t-)  very  well  for  (Th/m  g-1)  = 100, 

i.  £ O 

but  the  prediction  of  E(t2)  becomes  less  accurate  for  lower  acceleration  rates, 
and  is  over-estimated  by  a factor  of  2 at  (Th/m^g-1)  - 10.  However,  the  lift- 
off time  t2  remains  quite  accurate  over  the  considered  range.  It  Is  also  ap- 
parent from  this  figure  that  the  (inviscid)  energy  release  given  by  Equation 
(30)  can  be  represented  by  a power  law,  i.e..  Equation  (2),  with  6 “ 1 up  to 
a time  t^  such  that  (Th/m^g)  (t^/Isp)  <<  1 . A reasonable  approximation  is 


t < t,  ■ 0.1  Isp/(Th/m  g) 
— ( o 


6 ■ < (t2“*)/(t2’^P  “ ^2 


0 


t t2  • 1.72  Isp/(Th/m^g)  J 


(36) 
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FIGURE  7.  IMVISCIO  APPROXIMATION  FOR  ENERGY  DEPOSITED  INTO  THE  UNSTEADY 
PLUME  AT  GROUND  LEVEL 
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^ However,  these  results  must  now  be  qualified  by  examining  the  effect  of  vis- 

cosity, which  has  been  shown  to  predominate  the  mass  flow  results.  Equation 
^ (26)  has  been  evaluated  using  the  previously  cited  steady  plume  solution  with 

turbulent  mixing^  to  produce  the  curves  displayed  in  Figure  (8).  Comparison 
of  Figures  (8)  and  (6)  indicates  that  the  initial  growth  of  energy  flux,  E,  and 
Its  subsequent  rapid  cut-off  are  largely  the  result  of  the  corresponding  grov/th 
and  cut-off  of  the  mass  flow  rate.  These  curves  have  been  integrated  to  obtain 
1 the  total  energy  released  E(t),  corresponding  to  Equation  (30).  The  viscous 

and  inviscid  results  are  compared  in  Figure  (9).  It  is  clear  that  mixing  re- 
i duces  dramatically  both  the  total  amount  of  energy  deposited  at  the  launch  site, 

E(t2),  and  the  lift-off  time,  t^i  as  well  as  increasing  the  rate  of  energy 
deposition;  g > 1 for  t < t2* 
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The  viscous  calculations  were  not  carried  to  higher  accelerations,  i.e., 

Th/m^g  >>  2,  because  the  nozzle  exit  plane  would  be  closer  to  the  ground  at 
time  t2  than  the  first  axial  station  of  the  steady  plume  solution^  selected  for 
integration.  The  first  axial  station  selected  was  that  at  which  the  turbulent 
mixing  layer  just  reached  the  plume  axis;  thus  viscous  effects  would  diminish 
significantly  at  stations  closer  to  the  nozzle  exit  plane.  Therefore,  it  is 
clear  that  the  inviscid  solutions  must  represent  the  actual  behavior  for  suf- 
ficiently high  acceleration  rates,  e.g..  Th/m  g > 20.  The  dependence  of 
• o 

E/<E  > on  Th/m  g will  therefore  be  controlled  by  turbulent  mixing  for  low  ac- 
o o 

celeratlon  rates.  Th/m  g < 2,  but  will  be  essentially  inviscid  at  high  acc6ler- 

o . 

ation  rates.  Th/m  g > 20.  A minimum  in  the  ratio  E/<E„^  may  be  expected,  as 
o 

suggested  in  Figure  (10),  for  acceleration  rates  In  the  Intermediate  range 

2 < Th/m  g < 20.  However,  the  transition  curve  from  the  viscous  to  inviscid 
o 

solution  shown  in  this  figure  is  conjectured. 

3.5  Representative  Values  of  Rq  and  t2  “ To  demonstrate  some  repre- 

' e 

sentative  values,  consider  a specific  Impulse  Isp-250s  and  a thrust  Th  - 2x10  N.. 
As  a very  low  acceleration  lift-off  example,  take  Th/m^g  ■ 1.17.  for  which  the 
viscous  solution  clearly  applies.  Assume  p A » Th/10  and  y ■ 1.20,  giving  an 

• ® O Q 

initial  (nozzle  exit  plane)  energy  flux  Eq  ■ 4.9x10^  N*m/s  and  a total  energy 

deposition  at  the  ground  surface  E • 2x1o''  N*m.  Thus  R * 157m  and  t«  “ 12s  In 
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FIGURE  8 


ENERGY  DEPOSITION  RATE  INTO  THE  UNSTEADY  PLUME  AT  GROUND  LEVEL 
(INCLUDING  TURBULENT  MIXISS)  • 
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this  case.  On  the  otherhand,  as  a very  high  acceleration  lift-off  example, 
take  Th/m  g = 100,  for  which  the  inviscid  solution  should  apply.  Again  assume 

O • Q 

p A “ Th/10  and  y = 1.20,  giving  <E>  - A. 9x10  N-m/s  as  in  the  previous  case. 

C C Q 

However,  E = 8.8x10  N-m,  R = 56  m and  t.  = 5s  in  this  case. 

o 2 
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IV.  ASYMPTOTIC  I NV I SC  ID  STRUCTURE 


4.1  Unsteady  Mach  Disc  and  Plume  Collapse  - As  noted  by  Erdos  and 
Del  Guidice^  , the  supersonic  plume  is  quasi-steady  from  the  sonic  line  to 
the  Mach  disc  in  the  sense  that  it  can  be  accurately  described  at  any  instant 
by  a steady  solution  subject  to  the  boundary  conditions  prevailing  at  the 
sonic  line.  However,  the  motion  of  the  Mach  disc  is  controlled  by  the  instantan- 
eous pressure  imposed  on  its  downstream  side;  thus  its  position  is  determined 
by  the  interaction  of  the  quasi-steady  plume  with  the  unsteady  blast  field. 

This  interaction  in  general  requires  numerical  solution  of  the  unsteady  shock 
layer  contained  between  the  Mach  disc  and  blast  wave^ . However,  considerable 
insight  regarding  the  role  of  the  unsteady  shock  layer  can  be  gained -wi thout 
reverting  to  the  numerical  solution  itself.  For  example,  as  pointed  out  by 
Schmidt  and  Shear®,  at  late  times,  when  the  blast  wave  weakens  and  the  pres- 
sure acting  on  the  Mach  disc  decays  to  atmospheric,  the  plume  and  blast  fields 
uncouple.  Thus,  during  the  collapse  of  the  plume,  the  instantaneous  Mach  disc 
position  can  be  approximated  closely  by  the  steady  state  correlation 
formula 

i 

2^/d  = 0.70  (yp^/pj  M^  (37) 

where  d is  the  diameter  of  the  gun  barrel,  Pg(t)  is  the  muzzle  exit  pressure 
and  M^(t)  is  the  muzzle  exit  Mach  number  (usually  unity).  This  behavior  is 
demonstrated^  in  Figure  (11). 

In  the  case  of  a rocket  launch  this  formula  can  be  rewritten  as: 


r 

m 


1 .40 


^ P«c 


(33) 


where  c Is  the  exhaust  sound  speed  at  the  sonic  line  and  E is  given  by  Equa- 
tion (26).  This  form  clearly  shows  the  dependence  of  the  Mach  disc  positior, 
during  the  time  of  plume  collapse,  on  the  rate  of  energy  desposition,  rather 
than  the  amount. 
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On  the  other  hand,  during  the  early  period  when  the  plume  is  growing,  rather 
than  collapsing,  the  pressure  imposed  on  the  Mach  disc  is  essentially  that 
produced  by  the  blast  wave.  Thus,  a correlation  formula  of  the  form  given 
by  Equation  (10)  may  be  expected  to  apply;  this  is  confirmed^  by  the  numeri- 
cal examples  shown  in  Figure  (12).  However,  it  should  be  noted  that  the 
range  of  applicability  of  Equation  (10)  must  be  limited  to  the  duration  of 
strong  coupling  between  the  blast  field  and  plume.  Interestingly,  as  indi- 
cated in  Figure  (12),  the  correlation  appears  to  begin  to  deteriorate  when 
the  forward  motion  of  the  Mach  disc  drops  below  Mach  2 (based  on  atmospheric 
sound  Sfeed). 

Therefore,  the  behavior  of  the  Mach  disc  over  the  entire  duration  of  exist- 
ence of  the  unsteady  plume  can  be  approximated  by  application  of  the  corre- 
lation formulas,  or  scaling  laws,  represented  by  Figure  (12)  and  Equation 
(37)  or  (38),  in  their  appropriate  range  of  validity,  i.e.,  t <<  t2  and 
t -*■  tj.  respectively. 

4.2  Asymptotic  Blast  Field  - The  contact  surface  and  blast  wave 
trajectories  at  early  times,  i.e.,  t < t^,  have  been  shown^  in  Figures  (4) 
and  (5)  to  correlate  in  the  form  of  Equations  (8)  and  (9)*  During  this  time 
the  unsteady  shock  layer  flowfield  has  been  described  numerically^.  No 
attempt  to  correlate  the  radial  distribution  of  flow  properties  was  made^  be- 
cause of  the  ease  and  speed  with  which  the  numerical  solution  could  be  car- 
ried out.  However,  as  the  distance  between  the  shock  and  contact  grows, 
and  the  strength  of  the  blast  wave  diminishes,  the  accuracy  of  the  numerical 
solution  deteriorates,  and  recourse  to  an  asymptotic  formulation  becomes 
desirable^.  Note  that  for  t > tj,  the  plume  has  collapsed,  and  the  exhaust 
gases  remaining  at  the  launch  site  are  contained  in  a constant  (atmospheric) 
pressure  volume  bounded  by  the  inviscid  contact  surface.  However,  the  blast 
wave  continues  to  travel  outward  producing  a local  perturbation  in  pressure, 
temperature,  and  velocity  in  the  atmosphere.  The  asymptotic  solution  for 
the  flowfield  produced  by  such  a spherical  shock  wave  has  been  derived  by 
Whitham**,  Ranlet  and  Erdos^  have  correlated  their  numerical  solutions  for  a 

V 

particular  weapon,  the  M16  rifle,  at  late  times  when  the  blast  shock  Mach  | 

number  is  between  1.10  and  1.04,  In  terms  of  the  variables  from  Whitham's 

i ' 
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asymptotic  solution,  namely: 


i - 0 = f/n) 


(39) 


L 1)  = if^(n) 


m 


(J)  - 1 - --)  = f,(n) 

.d  p c 2 


(^♦1) 


OD  00 


\ 


I 


\ 

i 


1 

I 

i 

( 

I 

r 

c 

I 


where 


r - a t 

00 

d Jin  (r/dT 


('♦2) 


which  applies  to  ->  I and  -►  ”,  [as  opposed  to  the  strong  blast  theory 
used  previously  (Equations  2 - 10),  which  is  limited  to  » 1.]  The 
asymptotic  shock  solution  was  stated  earlier.  Equation  (11).  (Note  that  the 
gun  bore  diameter,  d,  was  used  by  Ranlet  and  Erdos^  as  the  length  scale, 
rather  than  R^.  However,  this  simply  introduces  a constant  factor  which  is 
Irrelevant  to  the  following  discussion).  It  should  be  pointed  out  that 
Whitham**  finds  the  asymptotic  flowfield  to  be  an  N wave,  i.e.,  a dual  shock 
system  with  the  functions  f^  and  f^  linear  and  quadratic,  respect ively,  be- 
tween the  two  shocks. 

The  pressure  distribution  across  the  shock  layer,  from  contact  to  blast  wave, 
obtained^  from  the  numerical  solution  for  an  Ml6  rifle.  Is  plotted  in  Figure 
(13)  in  the  form  given  by  Equation  (39).  It  Is  evident  that  the  numerical 
results  are  well  correlated  in  terms  of  the  variable  h,  although  a true  N 
wave  is  not  obtained.  The  corresponding  density  distributions  are  shown  in 
Figure  (1^).  In  this  figure  the'presence  of  an  entropy  layer  adjacent  to  the 
contact  surface  can  be  clearly  seen;  however,  outside  the  entropy  layer  the. 
density  distributions  are  also  well  correlated.  The  entropy  layer  consists 
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of  shock-heated  air  generated  at  early  times,  when  The  existence  of  such 

an  entropy  layer  is,  of  course,  not  included  in  Whitham's  asymptotic  theory. 

The  difference  between  the  density  and  velocity  in  the  form  given  by  Equation 
{k)  is  shown  in  Figure  (15).  In  this  figure  the  effect  of  the  entropy  layer 
appears  even  more  pronounced;  however,  note  that  the  considered  velocity- 
density  parameter  is  a higher  order  quantity  than  either  the  pressure  incre- 
ment or  density  increment.  It  should  also  be  noted  from  Figures  (l4)  and  (15) 
as  well  as  Figure  (13).  that  the  solutions  neat  the  blast  wave  correlate  at 
the  earliest  times  considered  and  the  correlated  region  subsequently  expands 
i nwa  rd . 

Since  the  asymptotic  behavior  of  the  blast  field  is  independent  of  the  initia- 
ting mechanism,  the  numerical  solutions  for  muzzle  blasts  of  widely  varying 
guns  should  also  correlate  (asymptotically)  in  terms  of  the  parameters  given 
by  Equations  (39)  " (^2).  This  has  been  verified  by  a comparison^  of  numerical 
solutions  for  the  M16  rifle  (5.56  mm  bore),  a special  purpose  infantry  weapon, 
SPIW  (a  5.8  mm  rifle  with  a non-standard  round)  and  a howitzer  (155  mm  bore). 

As  can  be  seen  in  Figures  (16)  and  (17)  the  pressure  and  density  distributions 
for  these  three  weapons  are  again  well  correlated.  However,  the  parameter 
representing  the  difference  between  the  density  and  velocity  (Equation  4l)  re- 
tains a dependence  on  the  initiating  mechanism  which  does  not  appear  to  vanish 
aymptotical ly.  However,  as  shown  in  Figure  (18),  scaling  of  this  parameter  to 

a reference  value  chosen  as  the  value  at  r = a t,  minimizes  the  dependence  on 

00 

the  Initiating  mechanism  and  produces  a satisfactory  correlation  for  the  three 
weapons. 

Although  these  correlations  were  originally  developed  to  demonstrate  that  the 
numerical  computation  did,  in  fact,  retain  a satisfactory  degree  of  accuracy 
as  ->  1 and  r^  •>  «>,  provided  adequate  grid  control  was  used^,  they  furnish 
additional  valuable  insight  regarding  the  asymptotic  behavior  of  far  field. 

This  is  particularly  useful  in  consideration  of  the  boundary  conditions  on  the 
turbulent  mixing  field  which  develops  in  the  presence  of  the  entropy  layers 
bounding  both  sides  of  the  inviscid  contact  surface. 
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V.  NONEQUILIBRIUM  UNSTEADY  SHOCK  LAYER 

5.1  Introductory  Remarks  - The  unsteady  shock  layer  is  defined  as 
the  region  between  the  plume  Mach  disc  and  the  blast  wave.  It  contains  pro- 
pellant gases  at  temperatures  on  the  order  of  combustion  chamber  temperature, 
located  between  the  Mach  disc  and  the  contact  surface,  and  shock-heated  air, 
located  between  the  contact  surface  and  blast  wave  (but  mainly  adjacent  to 
the  contact  surface).  Turbulent  mixing  of  the  propellant  gases  and  air  a- 
long  the  contact  surface  will  be  considered  in  a following  section;  for  the 
present  they  will  be  assumed  to  be  separated  by  the  contact  surface,  an  im- 
permeable membrane  in  the  context  of  an  inviscid  solution. 

Chemical  reactions  occurring  in  the  combustion  chamber  typically  freeze 
during  the  rapid  expansion  in  the  exhaust  nozzle,  and  remain  frozen  in  the 
supersonic  portion  of  the  plume  within  which  a continued  expansion  of  the 
flow  occurs.  Chemical  reactions  may  be  expected  to  reinitiate  in  the  stag- 
nation region  where  the  exhaust  flow  impinges  on  the  ground  surface  and  turns 
radially  outward.  However,  a second  freeze-out  is  likely  in  the  supersonic 
plume  which  forms  along  the  ground.  Figure  (2).  For  the  purposes  of  the 
present  exploratory  study,  the  impingement  region  has  been  neglected,  and 
all  chemical  reactions  have  been  assumed  to  freeze  in  the  exhaust  nozzle  and 
remain  frozen  until  the  unsteady  shock  layer  is  entered.  (However,  this  as- 
sumption will  be  subsequently  reexamined  in  the  light  of  computational  re- 
sults.) 

Supersonic  expansion  of  the  propellant  gases  produces  the  well-known  under- 
expanded plume  configuration  shown  in  Figure  (1),  including  the  barrel  shock 
needed  to  match  the  continuously  decreasing  pressure  within  the  core  of  the 
plume  to  the  ambient  field  prevailing  along  its  lateral  boundaries.  Conse- 
quently, the  supersonic  expansion  is  not  a true  source  flow,  and  the  differ- 
ences, even  along  an  axis  of  symmetry,  can  produce  an  appreciable  shift  in 
the  Mach  disc  location,  as  demonstrated  in  References  (1)  and  (2).  Neverthe- 
less, the  supersonic  source  flow  model  is  extremely  useful  in  an  exploratory 
study  since  the  solution  is  provided  by  well-known  isentropic  flow  relations 
rather  than  by  a cumbersome  numerical  solution.  However,  Its  limitations  in 


accuracy  ard  its  lack  of  dependence  on  the  ambient  pressure  should  be  recogniz- 
ed . 


The  flow  pattern  indicated  in  Figure  (l)  may  begin  at  earliest  times  with 
roughly  cylindrical  symmetry  (relative  to  the  launch  axis);  however,  a more 
nearly  spherical  symmetry  should  be  expected  near  the  blast  wave  after  a short 
time  (during  which  the  vehicle  does  not  move  appreciably).  The  success  with 
which  the  gun  blast  can  be  modelled  as  a spherically  symmetric  shock  layer  ^ 
suggests  the  same  approximation  should  be  equally  useful  in  this  case.  There- 
fore, the  supersonic  plume  is  represented  by  a spherical,  quasi-steady  source 
flow.  The  unsteady  subsonic  flow  between  the  Mach  disc  and  the  blast  is  des- 
cribed by  a finite-difference  solution  of  the  one-dimensional  time-dependent, 
inviscid  equations  of  motion  (with  spherical  symmetry),  in  the  same  manner  as 
employed  in  Reference  (l)  and  (2).  The  resulting  solution  should  be  most  ac- 
curate along  the  ground  (or  sea)  surface,  and  least  accurate  along  the  launch 
axis  (where  the  blast  is  considerably  less  intense).  In  addition,  the  existence 
of  turbulent  mixing  of  the  propellant  gas  and  air  along  their  interface  and  of 
buoyant  forces  is  pointed  out;  however,  neither  ii  included  in  the  present  in- 
viscid model.  The  turbulent  mixing  problem  will  be  addressed  in  the  following 
section. 


5.2  Governing  Equations  and  Solution  Algorithm  - The  system  of  govern- 
ing equations  (for  a perfect  gas)  and  the  solution  algorithm  have  been  described 
in  Reference  (l)  and  (2),  therefore,  they  will  be  only  briefly  reviewed  here. 

The  principal  modifications  relate  to  the  non-equil  ibrlum  chemical  kinetics, 
which  will  be  described  in  depth.  The  system  of  governing  equations  is  stated 
as: 


Continui ty: 

If  + ^ (pu)  - - k 

3t  3r  ^ r 

Momentum: 

^ (pu^p)  - - k 

3t  3r  ^ ^ r 

Energy: 

De  , Op  ' - dr 

5T*  " • »"  dT' 

m 

m 

(^*5) 
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Species :* 

State: 


Da . 0) . , 

dT  = T = ).  on  u 


p = pR  T Z a./M. 

O II 


(46) 

(47) 


e = h - - Z a.h  - R T Z a./M. 

p I ! o II 


h , = a.  + b.  T + c.  T 

ill  I 


(48) 


where  k is  an  integer  which  is  0 for  plane  symmetry,  1 for  cylindrical  sym- 
metry and  2 for  spherical  symmetry.  Only  k = 2 is  used  in  the  present  study. 
The  region  between  the  Mach  disc  and  the  blast  wave  is  divided  into  two  layers; 
the  first  extends  from  the  Mach  disc  to  the  contact  surface,  which  as  previously 
mentioned,  is  an  impermeable  interface  between  the  propellant  gases  and  air, 
and  the  second  extends  from  the  contact  surface  to  the  blast  wave.  These 
three  boundaries  are  treated  as  moving  surfaces  of  d i scont inui ty  across  which 
certain  gas-dynamic  jumps  occur.  Denoting  their  instantaneous  positions  as 
r^ , r^  and  r^,  respectively,  a stretched  coordinate  Cj  for  each  layer  is  de- 
fined by: 


i 


(r 


•ri)/(ri^^-ri) 


(50) 


where  1 * 1 refers  to  the  first  layer  (from  Mach  disc  to  contact)  and  i = 2 
to  the  second  (from  contact  to  blast  wave).  The  above  system  of  equations 
Is  transfromed  from  (r,t)  to  (6j»t)  coordinates  by: 


^ _ J_  _3_  (51) 

3r  “ Atj  3^1 

*Note  that  the  use  of  Equation  (46)  for  al 1 species  together  with  Equation  (43) 
produces  a redundant  system.  The  redundancy  is  used  to  obtain  an  independent 
check  on  error  accumulation.  The  error  defined  as 

I 

|1  - Z a I 
-5  ’’’ 

is  less  than  10  In  all  results  obtained  thus  far. 
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_3_ 

at 


_a_ 

at 


(l-C.)  u,  . 


C,  u 


i + 1 


Ar . 


35: 


where 


Ar.  = r. 


i + 1 


u. 

I 


dr. 



dt 


(52) 


(53) 

(5i|) 


The  two  domains  are  each  divided  into  a series  of  equally  spaced  intervals 

of  size  A^..  The  magnitude  of  A^^  and  are  independently  controlled  such 

that  the  magnitude  of  the  grid  size  (Ar.  ACj)  remains  within  a preselected 

range.  This  is  accomplished  by  adding  or  subtracting  a single  grid  point 

when  the  size  falls  out  of  the  range,  shifting  the  other  points,  and  linearly 

interpolating  the  data  at  the  new  grid  noint  locations.  Since  the  resulting 

2 

shift  in  grid  location  is  of  the  order  AC.,  no  loss  of  accuracy  is  incurred 
by  the  linear  interpolation. 


Equations  (43),  (44)  and  (46)  are  expressed  in  finite  difference  form  and 
solved  at  each  of  the  interior  grid  points  by  the  MacCormack  algorithm^.  Equa- 
tion (45)  is  solved  by  a Lagrangian  ;nethod,  i.e.,  the  streampath  dr/dt  = u is 
traced  from  a grid  pont  back  to  its  intersection  with  the  previous  time  plane. 
Recall  that  in  the  case  of  a perfect  gas  Equation  (45)  reduces  to: 


DS 

Dt 


0,  or  S 


constant , 


u 


(55) 


Thus,  interpolation  of  the  value  of  the  entropy, S,  from  the  point  of  inter- 
section of  the  streampath  with  the  previous  time  plane  provides  the  "exact" 
solution  in  this  case,  and  accordingly  any  numerical  diffusion  of  an  entropy 
gradient  is  minimized.  To  maintain  a consistent  degree  of  accuracy  in  the 
solution  of  energy  equation  with  chemical  reactions.  Equation  (45)  was  solved 
in  the  same  fashion  as  Equation  (55)* 
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e + p/p  = constant  on  — = u (5^ 

at 

where  p is  the  average  pressure  along  the  considered  streampath  segment. 

Solution  of  Equations  (^3)  and  by  the  MacCormack  algorithm  is  stable  (with- 

out any  additional  numerical  damping)  if  the  standard  CFL  cr i ter  ion  i s observed^ . 
In  simple  terms,  the  maximum  permissible  time  step  is  limited  by  the  time  for 
a wave  to  reach  an  adjacent  grid  point.  However,  addition  of  Equation  to 

the  system  introduces  a second  time  scale,  namely  that  for  chemical  relaxation. 

By  analogy  with  the  results  obtained  by  Cheng^°  for  combined  wave  and  diffusion 
processes,  the  time  step  criterion  for  the  present  system  is  postulated  as: 


speed 


where  the  character i St i c time  scales  are: 


chemical 
rel axat ion 


(u  + C|:)/Ar 


speed 


chemical 

relaxation 


a “i 
— (— ) 
3a 


and  the  frozen  sound  speed  c^  is  defined  by 


(ill) /(it  . i) 

'3p^''3p  p' 


which  can  be  analytically  evaluated  using  the  state  relations  provided  by 
Equations  (1*7),  (i*8)  and  (i*9) . 
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It  should  be  pointed  out  that  only  in  the  case  where  At  and  At  , . , 

wave  chemical 

are  of  distinctly  different  orders  of  magnitude  does  Equation  (57)  corre- 
spond to  using  the  minimum  of  the  two  time  scales.  On  the  other  hand,  if 
they  are  approximately  equal.  Equation  (57)  yields: 

At  . . , , or  At  = At 

chemical  2 wave  (61) 

relaxation  speed 


This  behavior  may  account  for  the  empirical  proportionality  constant  which 
Anderson^^  found  to  be  necessary  when  using  the  minimum  of  the  two  time 
scales.  The  case  in  which  the  time  scales  are  of  the  same  order  of  magni- 
tude corresponds  to  a true  non-equilibrium  situation;  the  cases  in  which 
they  are  distinctly  different  correspond  to  frozen  and  equilibrium  limits. 


At  the  bounding  surfaces  of  discontinuity,  namely;  the  blast  wave,  contact 
surface  and  Mach  disc,  the  method-of-character i st ics  procedure^  for  determi- 
nating the  motion  of  these  surfaces  and  flow  properties  thereon  has  been  re- 
tained using  a locally  defined  frozen  sound  speed  and  isentropic  exponent. 

The  effects  of  chemical  reaction  are,  therefore,  neglected  within  narrow 
zones  of  thickness  (u  ^ c^)  At  adjacent  to  each  of  these  surfaces.  The  isen- 
tropic exponent,  y , is  evaluated  locally  from: 

y = '’4  p/p  (62) 


As  indicated  in  the  introductory  remarks,  the  flowfield  upstream  of  the  Mach 
disc  is  approximated  by  a supersonic  spherical  source.  In  effect,  the  sonic 
line  (surface)  is  located  at  a radius  r*,  and  the  Mach  number  at  any  radius 
r > r*  is  obtained  from  the  isentropic  flow  relationship. 


(r*/r)^ 


M [1  + ^ 


Y+1 

2(y-1) 


(63) 


The  instantaneous  conditions  on  the  upstream  side  of  the  Mach  disc  are  ob 
tained  from  the  above  equation  with  r ■ r^  (t)  and  with  the  entropy,  total 
enthalpy  and  species  concentrations  on  a streampath  traced  from  the  sonic 
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line  to  a position  rj(t).  The  latter  procedure  introduces  a tir^e  lag  between 
the  variation  of  sonic  line  conditions  and  the  resulting  variation  of  condi- 
tions at  the  Mach  disc  which  is  only  appreciable  when  rj  >>  r*. 


Ambient  atmospheric  conditions  are  assumed  to  prevail  upstream  of  the  blast 
wave  (i.e.,  any  effects  of  precursor  waves  due  to  ejection  from  a silo,  for 
example,  are  neglected).  The  instantaneous  velocities  of  the  Mach  disc  and 
blast  wave,  u^  and  u^,  are  determined  by  matching  the  solution  of  the  Rankine- 
Hugoniot  shock  Jump  relations  to  the  characteristic  compatibility  equation  a- 
long  the  wave  impinging  each  shock  from  the  interior  (i.e.,  an  upstream  run- 
ning wave  at  the  Mach  disc  and  a downstream  running  wave  at  the  blast  wave). 
The  velocity  of  the  contact  surface,  u^,  is  determined  by  equating  the  pres- 
sures and  (normal  component  of)  velocities  at  the  contact  given  by  the  com- 
patibility equations  on  the  pair  of  waves  impinging  on  the  contact  (one  from 
the  propellant  side  and  one  from  the  air  side).  Due  to  the  assumption  of 
spherical  symmetry,  no  tangential  component  of  velocity  exists.  Therefore, 
only  the  entropy  or  related  thermodynamic  properties,  and  the  species  con- 
centrations are  discontinuous  across  the  contact  surface. 


5.3  Gas  Dynamic  Properties  - As  a prelude  to  analysis  of  the  chemi- 
cal kinetics  of  the  propellant  gases  in  the  blast  field,  it  is  instructive  to 
study  the  gas  dynamics  of  the  problem  with  a view  toward  characterizing  the 
time  scales,  as  well  as  temperature,  pressure  and  velocities.  This  prelimi- 
nary assessment  can  be  conveniently  carried  out  on  the  basis  of  a binary  per- 

2 Z o 

feet  gas  model;  one  gas  represents  air  (y  “ l.A  and  R ■=  1716  ft^/sec^/  R)  and 

2 Z o 

the  other  represents  the  propellant  gas  (y  » 1.25  and  R ■ 211*0  ft  /sec  /°R) . 

In  this  case  Equation  (55)  is  used  in  place  of  Equation  (56),  iLj  “ 0,  etc. 

As  pointed  out  earlier,  the  energy  deposition  process  can  be  veiwed  as  be- 
ginning impulsively,  t.  <<  tj^,  and  the  mass  and  energy  flux  into  the  blast 
field  can  be  approximated  by  integrating  flow  properties  radially  across  a 
steady  exhaust  plume  at  various  stations  behind  the  nozzle  exit  plane.  It 
is  assumed  that  a quasi-steady  supersonic  exhaust  flow  and  plume  will  be 
established  behind  the  leading  (double)  shock  system  expelled  from  the  nozzle 
as  the  chamber  pressure  rises;  thus  the  time  to  establish  a quasi-steady  plume 
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should  be  even  less  than  t..  (Shadowgraphs®  of  the  muzzle  blast  of  the  Ml6 
rifle  when  the  projectile  is  still  within  one  bore  diameter  of  the  muzzle 
corroborate  this  assumption.)  A quantitative  assessm^t  of  the  effect  of  the 
relative  magnitude  of  t.  on  the  blast  field  processes  of  present  interest  has 
been  carried  out  by  comparing  the  solution  for  an  "instantaneous"  ignition 
(t.  = O)  with  solutions  for  a finite  ignition  time.  The  ignition  process  has 
been  represented  by  a chamber  pressure  rise  time  curve  taken  from  Reference 
(5).  Using  the  engine  and  propellant  parameters  selected  for  Case  1 in  Ref- 
erence (5),  e.g.,  axial  engine  length  of  10  ft.,  nozzle  area  ratio  of  ^.0, 

etc.,  the  steady  chamber  pressure  is  reached  in  about  25  milliseconds.  The 
pressure  at  the  sonic  line  was  assumed  to  respond  instantaneously  to  the 
chamber  pressure  history.  The  sonic  line  pressure-t ime  curve  is  reproduced 
here  as  Figure  (19).  An  intermediate  case  in  which  t^  = 2.5  milli- 
seconds has  also  been  carried  out,  using  the  same  form  of  the  pressure  rise 
curve. 

The  gas  temperatures  on  each  side  of  the  Mach  disc  and  at  the  blast  wave  are 
shown  in  Figure  (20)  as  a function  of  time.  The  temperatures  on  each  side 
of  the  contact  surface  are  shown  in  Figure  (21).  It  is  evident  that  sub- 
stantial differences  in  the  shock  layer  temperature  field  from  the  Mach  disc 
to  the  contact  occur  during  the  ignition  transient,  i .e. , t < tj,  but  for 
t > tj  the  solutions  tend  to  converge.  This  results  becomes  obvious  as  the 
time  for  energy  deposition  becomes  small  compared  to  the  observation  time. 

Thus  the  effect  of  finite  ignition  time  need  be  considered  only  in  (extreme) 
cases  where  it  is  comparable  to  the  time  scale  of  interest. 

Radial  distributions  of  temperature  at  various  instants  are  displayed  in  Fig- 
ure (22).  It  Is  clear  that  the  bulk  of  the  propellant  gas  is  at  a tempera- 
ture characterized  by  that  behind  the  Mach  disc,  which  in  this  perfedt  gas 
example  is  close  to  8000°R.  The  blast  wave  produces  a very  minor  heating  of 
the  air,  although  the  volume  of  air  processed  by  the  blast  wave  grows  quite 
rapidly.  The  temperature  ratio  across  the  contact  surface  is  of  the  order  of 
6 to  1 , and,  as  will  be  discussed  later,  may  produce  a thermally-driven  tur- 
bulent mixing  along  this  interface.  A steep  temperature  gradient  is  also  ob- 
served on  the  propellant  side  of  the  contact,  which  may  be  expected  to  further 
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FIGURE  20.  THE  TEMPERATURE  AT  THE  MACH  DISC  AND  BLAST  WAVE  AS  A FUNCTION  OF  TIME 
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contribute  to  the  turbulent  mixing  of  hot  propellant  gases  and  air. 

Radial  distribution  of  gas  pressure,  temperature  and  veiocity  at  a much 
later  time  are  shown  in  Figures  (23),  (2^)  and  (25).  The  assumption  of  a 
chemical  freeze  in  the  rocket  nozzle  and  supersonic  exhaust  plume  are  justi- 
fied heur i St ical ly  on  the  basis  of  the  rapid  pressure  and  temperature  drop 
and  increasing  gas  velocity  which  characterize  the  exhaust  nozzle  and  super- 
sonic plume,  cf.  Figures  (23),  (24)  and  (25).  These  considerations  pertain 
to  the  quasi -steady  flow  which  exists  upstream  of  the  Mach  disc.  However, 
the  present  unsteady  analysis  is  initiated  when  the  Mach  disc  first  forms, 
which  is  at  a position  close  to  the  sonic  line  (as  shown  in  Figure  l).  At 
this  early  time,  the  assumption  of  a chemical  freeze  upstream  of  the  Mach 
disc  is  less  tenable  and  the  initial  conditions  may  more  closely  correspond 
to  the  local  equilibrium,  as  will  be  seen  below. 

The  nearly  constant  pressure  plateau  between  the  Mach  disc  and  contact  sur- 
face exhibited  in  Figure  (23)  extends  well  into  the  region  between  the  contact 
and  blast  wave  at  later  time;  cf.  Figures  (13)  and  (16).  Thus,  the  velocity 
field,  shown  in  Figure  (25),  is  correspondingly  constant.  The  fact  that  no 
velocity  jump*  exists  across  the  contact  surface  is  reemphasized,  in  contrast 
to  the  temperature  jump,  seen  in  Figure  (24),  which  persists  forever  in  the 
absence  of  turbulent  motion.  Therefore,  the  asymptotic  behavior  of  the  pro- 
pellant gas  cloud  must  clearly  be  dominated  by  thermal  rather  than  inertial 
or  pressure  gradients. 

5.4  Thermo-chemical  Properties  - The  following  chemical  reaction 
system  was  selected  as  a representative  mechanism  for  afterburning  of  a hypo- 
thetical hydrocarbon  fuel  which  produces  a combination  of  H^O,  H21  CO2  and  CO 
and  N2  as  principal  exhaust  products. 


*The  assumption  of  spherical  symmetry  precludes  the  existence  of  a finite  velo 
city  component  parallel  to  the  contact  surface,  which  could  have  a jump. 
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Rates  for  the  collisional  dissociation  of  H2O  and  O2  by  N2  impact  were  taken 
from  Reference  (12);  the  remaining  rates  were  taken  from  Reference  (13).  The 
thermodynamic  data  for  the  const ! tuent  species  was  taken  from  Reference  (14) 
and  Reference  (15). 


For  il'ustrative  purposes,  the  chemical  composition  of  the  exhaust  stream  at 
the  nozzle  exit  plane  was  taken  from  the  study  by  Woodroffe^®.  Two  cases 
were  considered,  designated  as  fuel  rich  and  fuel  lean;  the  mole  fractions 
are  listed  in  Table  I.  As  mentioned  earlier,  the  composition  was  assumed 
to  be  frozen  from  the  nozzle  exit  plane  to  the  Mach  disc  of  the  ground  sur- 
face plume. 


TABLE  I 

Assumed  Chemical  Composition  of  Nozzle  Exhaust  Stream  By 
Mole  Fraction  (From  Woodroffe,  Reference  (I6) 

Species  Fuel  Lean  Fuel  Rich 


H20 

.30 

.10 

”2 

.05 

.25 

^^2 

.30 

.10 

CO 

.05 

.25 

"2 

.30 

.30 
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The  species  mass  fractions  at  the  propellant  side  of  the  contact  surface 
are  shown  as  a function  of  time  In  Figures  (26)  and  (27).  The  rapid  varia- 
tion during  the  first  10  psec,  together  with  the  fact  that  At  >>  At  ^ , 

chem  wave 

indicate  that  a near-equ i 1 i br i urn  cond  i t ion  preva i 1 s initially  and  the  antici- 
pated freeze-out  will  develop  later.  Therefore,  non-equilibrium  chemistry 
should  be  included  in  the  rocket  nozzle  and  supersonic  plume  solutions.  In 
addition,  it  may  be  more  appropriate  to  initialize  the  shock  layer  with 
equilibrium,  rather  than  frozen,  species  concentrations. 

The  variation  of  species  concentration  across  the  layer  from  the  Mach  disc 
to  the  contact  surface  at  t = 25^  psec  (for  the  fuel  lean  case)  is  shown  in 
Figure  (28).  The  corresponding  temperature  distribution  is  depicted  in  Fig- 
ure (29).  It  is  apparent  that  the  reaction  zone  is  just  downstream  of  the 
Mach  disc,  with  the  flow  equ i 1 i br i brat i ng  as  it  approaches  the  contact  sur- 
face. On  the  other  hand,  examining  the  complete  temperature  distribution 
across  the  shock  layer,  in  Figure  (30) , it  also  is  apparent  that  the  contact 
surface  is  artificially  separating  the  equilibrated  propellant  species  from 
a surrounding  "bath"  of  shock-heated  air.  Therefore,  additional  afterburning 
must  be  anticipated  along  the  contact  at  a rate  limited  only  by  turbulent 
dilution  of  the  mixture.  The  mixing  problem  is  addressed  in  the  following 
section.  The  present  non-equilibrium  chemical  kinetic  model  is  also  an 
idealization  in  the  sense  that  (a)  the  propellant  exhaust  products  are  as- 
sumed to  consist  of  a mixture  of  relatively  simple  molecules  such  as  water, 
hydrogen,  carbon  dioxide,  carbon  monoxide  and  nitrogen,  and  (b)  the  reaction 
scheme  only  includes  collisional  dissociation  of  water  vapor  and  oxygen  mole- 
cules by  molecular  nitrogen,  a few  representative  hydroxyl  radical  reactions, 
and  the  Zeldovich  reactions  for  atomic  nitrogen  production.  Future  effort 
should  focus  on  expansion  of  the  present  basic  reaction  scheme  to  include  all 
well-established  reactions  haveing  a potentially  measurable  impact  on  the  tem- 
perature and  concentration  of  optically  significant  species,  with  a view  toward 
quantifying  the  intensity  and  duration  of  the  optical  flash. 
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FIGURE  30.  TEMPERATURE  DISTRIBUTION  ACROSS  THE  ENTIRE  UNSTEADY  SHOCK  LAYER 
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VI . UNSTEADY  TURBULENT  MIXING  LAYER 

6.1  General  Approach  - Observations®  of  muzzle  blasts  show  the 
presence  of  intense  turbulence  at  very  early  stages  of  development  of  the 
flowfield,  which  eventualiy  engu I fs  the  ent i re  propellant  cloud.  Press  pho- 
tographs of  unclassified  missile  launches  also  exhibit  the  same  sort  of  tur- 
bilent  propel iant  cloud.  (In  retrospect,  it  could  have  been  very  interesting 
to  examine  movies  of  rocket  iaunches  to  study  the  details  at  very  early  times 
following  ignition.)  However,  the  source  of  the  turbulence  is  not  identifi- 
able from  these  observations.  Among  the  possible  sources  are  high  frequency 
pressure  osciliations  originating  in  the  combustion  chamber, shear-induced 
turbulence  along  the  lateral  boundaries  of  the  exhaust  piume  and/or  ground 
surface  plume,  and  Taylor  instability^^  of  the  contact  surface.  The  last 
mechanism  is  produced  by  acceleration  of  a light  gas  toward  a heavy  gas^^. 

Thus,  a spherical  volume  of  hot,  light  gas  (propellant  gases)  expanding  into 
a heavier  gas  (air)  should  experience  Taylor  instability  along  its  front  (the 
contact  surface).  As  previously  noted,  no  velocity  jump  exists  across  the 
contact,  but  an  extreme  temperature  Jump  is  present.  Therefore,  the  somewhat 
unconventional  notion  of  thermally-driven  turbulence  production  along  the 
contact  surface  is  suggested.  If  the  ideas  of  similitude  between  heat  and 
mass  transfer,  and  of  a turbulent  Prandtl  number,  representing  the  ratio  of 
eddy  viscosity  to  eddy  conductivity,  near  unity,  are  acceptable,  then  the  idea 
of  a similitude  between  an  eddy  viscosity  produced  by  velocity  gradients  and 
an  eddy  conductivity  produced  by  temperature  gradients  should  be  palatable. 
Although  temperature  fluctuations  are  conventional ly  regarded  as  the  product 
of  velocity  fluctuations,  it  is  suggested  herein  that  the  opposite  may  occur 
in  the  absence  of  shear  forces  as  a driving  mechanism  for  the  production  of 
velocity  fluctuations. 

The  development  of  turbulence  models  and  closure  statements  for  steady  bound- 
ary layer  and  wake  (free  mixing)  type  flows  has  been  the  subject  of  consider- 
able research,  cf.,  the  review  by  Mel  lor  and  Herring  (Reference  18).  While 
extension  of  these  models  to  unsteady  flows  is  not  conceptually  difficult  the 
experimental  data  base  is  considerably  more  meager.  Unfortunately,  most  at- 
tention in  the  study  of  unsteady  boundary  layers  has  focussed  on  definition 


-62- 


TR  2k6 


I 


i 

t 

f 

!' 


1 


f 


( 

I’'- 


[' 


and  prediction  of  unsteady  separation.  In  this  regard,  the  subject  problem 
provides  an  attractive  framework  for  the  study  of  free  turbulent  mixing  since 
it  is  (at  least  approximately)  two-dimensional  (i.e.,  time  and  radius)  and 
free  of  such  complexities  as  laminar  sublayers  and  separation.  Furthermore, 
as  shown  in  the  preceding  section,  self-similar  inviscid  solutions  describe 
the  outer  inviscid  layer  at  early  time  (strong  blast  theory)  and  at  late  time 
(Whitham's  asymptotic  theory),  which  may  be  comb ined  us  ing  the  method  of  match- 
ed asymptotic  expansions.  The  inner  inviscid  layer  (between  the  Mach  disc  and 
contact)  still  appears  to  require  numerical  solution,  although  5 posteriori 
examination  of  the  available  solutions  suggests  certain  simplifying  approxi- 
mations may  provide  an  adequate  representation.  Thus,  the  boundary  conditions 
on  the  mixing  layer  along  the  contact  surface  are  well  defined  both  in  terms 
of  the  present  numerical  solution  and  asymptotic  theories. 

Based  on  the  observations  of  Baronti  and  Miller  (Reference  19)  for  steady 
boundary  layers  and  of  Telionis  (Reference  20)  for  unsteady  boundary  layers, 
integral  methods  continue  to  provide  as  good  renditions  of  turbulent  flows  as 
do  the  more  numerically  detailed  finite  difference  solutions.  Therefore,  a 
two  layer  (i.e.,  two-strip)  integral  method  analysis  of  the  subject  free  mix- 
ing problem  is  strongly  suggested  by  both  physical  and  pragmatic  computational 
considerations.  The  physical  interpretation  of  the  proposed  two  layer  model 
Is  indicated  in  Figure  (31).  The  first  layer  6^  grows  radially  outward  from 
the  contact  surface  while  the  second  layer  62  grows  radially  inward.  The 
variations  of  elemental  fractions  of  "propellant"  (i.e.,  the  sum  of  all  pro- 
pellant species)  and  of  "air"  (i.e.,  the  sum  of  N2  and  O2  from  the  air)  are 
indicated  in  this  sketch.  The  temperature  distribution  should  follow  roughly 
the  curve  shown  for  the  "propellant"  mass  fraction.  It  is  emphasized  that  the 
mixing  layer  model  is  formulated  in  a coordinate  system  fixed  at  the  moving 
contact  surface;  thus  the  mean  velocity  is  zero  everywhere. 


As  mentioned  above,  the  available  data  base  for  mixing  layers  of  this  type  is 
rather  meager.  However,  comparison  shown  in  Figure  (32)  of  the  experimental 
observations  by  Schmidt  and  Shear®  of  the  growth  of  the  turbulent  front  along 
the  axis  of  a muzzle  blast,  with  the  inviscid  contact  surface  prediction  by 
Erdos  and  Del  Guidice^  provides  one  set  of  data  which  can  be  interpretted  as  the 
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growth  of  6^.  A corresponding  observation  of  62  is,  however,  totally  ob- 
scured by  the  propellant  cloud. 


6.2  rormulat ion  of  a Turbulence  Model 
energy  (TKE)  equation  forms  the  basis  of  most  of 
order  closure  models  now  in  use.  The  well-known 
of  this  equation  is  included  below: 


- The  turbulent  kinetic 
the  widely  accepted  higher- 
term-by-term  interpretation 


Dk 
^ Dt 

(convection) 


3U 

pu  V ^ 
(production) 


_3_ 
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(pv^k  + p'v'  ) 
(diffusion) 


(dissipation) 


where : 


1 (u'2. 


(65) 


(66) 


(The  need  for  a mean  velocity  gradient  U/  y to  drive  the  production  of  tur 
bulent  kinetic  energy  is  clearly  seen.) 


Spaulding's^^  two-equation  model  of  the  TKE  equation  is  one  of  t!ie  more  widely 
employed  closures: 
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where  u_  ■ C pk  /e  (69) 

T o 

The  k constants  Cq,  C2  and  Cj  have  been  selected  by  Spaulding^^  to  pro- 
vide the  best  agreement  with  the  data  base  for  steady  turbulent  boundary  layers 
(mainly  incompressible): 
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Obviously,  since  3U/3y  = 0 in  the  present  problem,  the  necessary  production 
mechanism  is  absent  in  the  TKE  equation,  and  therefore  also  in  the  Spaulding 
closure  model.  However,  application  to  the  thermal  energy  (or  heat  balance) 
equation  of  procedures  analogous  to  those  used  to  derive  the  TKE  equation 
yields^^’^^  a turbulent  thermal  energy  (TTE)  equation,  having  a similar  term- 
by-term  interpretation: 


2i 

Dt 

(convection) 


2p  v'e 

3y 

(production) 


(d i f fus ion) 


^ ^ ^ 2 (|i)  ) = 0 

3y 

(molecular  transport  and  smearing) 


where : 


(T  - T)' 


A two-equation  closure  model  can  be  defined  by  analogy  with  the  above  TKE 
closure  model : 


• 2 2 
26-  L fill 

P Dt  S’  3y' 


-)  + P(l( 


'l7>  - s ^ <4  1^'  * S-  ^ 


where: 


ky  - Cq  P C id  ) ‘ /^ 


Since  no  data  base  is  available  with  which  to  calibrate  the  constants  Cq,  , 

and  Cj,  Spalding's^^  values,  given  by  Equation  (70),  for  the  TKE  closure 
have  been  adopted.  If  the  assumed  analogy  between  turbulent  transport  of  heat 
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and  mass  exists,  these  values  should  be  at  least  approximately  correct. 


6.3  Method  of  Solution  - The  system  of  equations  consists  of  the 
statements  of  conservation  of  mass,  energy  and  species,  the  two-equation 
TIE  model,  and  an  equation  of  state.  For  the  present  discussion,  the  species 
equation  can  be  dropped  and  a single-component  perfect  gas  assumed  for  illu- 
strative purposes.  The  system  can  then  be  stated  consisely  by: 


L (f) 


where  L - ^ + v ^ - A (k^ 


f * [p.  T,  9 , 4i] 


A - [0,  p’\  p'\  C^p"^] 


(76) 

(77) 

(78) 

(79) 


B 


9v  1 9p 
3y’  pCp  9r’ 


(li  (>!)' 


- <^).  (c 


(11) 

P02  ^3y^ 


Cj(^^/02  )] 


(30) 


(Note  that  y is  measured  outward  from  the  contact  surface  and  v is  the  mean 
velocity  component  in  they  di  rect  ion  produced  by  the  turbulent  heat  conduc- 
tion. Since  spherical  symmetry  has  been  assumed  there  are  no  gradients  along 
the  contact  surface  or  velocity  component  parallel  to  it.) 


Following  standard  procedures,  the  first  two  components  of  the  system  repi 
sented  by  Equations  (76)  through  (80),  namely  the  continuity  and  energy  equa- 
tion, can  be  integrated  from  y ■ 0 to  y 
yielding; 


6^  (t)  and  from  y 


0 to  6^  (t). 


d5 


I 

dt 


^ d I M 

- — ; 0 fy 


(81) 
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and 


d6.  6.  . 

p.  T.  — — + ^ 


i i dt 


Cp  dt 


p.v.T.  + f- 
III  dt 


PT  dy  + 1 |I 


(32) 


where  ( ) denotes  y » 0 (i.e.,  at  the  contact  surface),  i'  = 1 and  2 denotes 
o 

values  at  y • 6^  and  y = 6^  respectively,  and  (3/3y)^  - (d/dy)^  “ 0 is  as- 
sumed, in  accord  with  standard  practice  for  first  order  boundary  layer  theory. 
Substituting  Equations  (8l)  into  (82)  gives: 


-i 


TT  / pT  dy 


5. 


_d 

dt 


pdy 


_L  ^ - 

Cp  dt 


H 

‘t  3y 
\ / 


(83) 


Utilizing  the  perfect  gas  equation  of  state  the  first  term  in  Equation  (83) 
becomes : 


dt 


p 

/ pT  dy  » 

dt 

Jo 

(84) 


Assume : 


P - P, 
^2  ■ 


y (1  + cos  irn) 


(85) 


where 


n ■ (y  - dj  ■ ^2^ 


(86) 


Then: 


d5 

dp 

L - 

/ \ 
fp,  +P2  1 

IP  dt  ^ ^i 

dti 

1 

^i  dt 

^ / 

T 

'l  dt  \ 2 


■ *2 


(87) 


' ^ . 3T 

“ Cp  dt  ■ 3y 
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If  p,  and  are  prescribed  functions  of  time,  from  the  inviscid  soiution, 
then  Equation  (87)  provides  two  equations  (i  = 1 and  2)  for  two  unknowns,  6^ 
and  62,  in  terms  of  the  turbulent  heat  conduction  rate  across  the  contact  sur- 
face, (k^  3T/3y)^.  Thus,  it  is  only  necessary  to  evaluate  k^  at  y = 0,  using 
the  last  two  components  of  the  system  represented  by  Equations  (76)  through 
(80).  However,  it  is  necessary  to  assume  distribution  functions  for  ¥2"  and 
i}i;  the  following  have  been  selected. 


(1  - cos  2n  n)/2 


(1  - cos  2n  n)/2 


Therefore,  the  two-equation  TTE  model  becomes: 


de^  4 2Tr^ 

o.  o o_  o 

dt  ■ p , ^2  p 3y 


- - C 


S’  2ir^  <1.^ 
o o 

’ ‘'o 


• !o_  /llf  c 6^/e  ^ 

'2  Po  ^ I'y/o  3 V ° 


where 


( — ) 

' 3y  o 


A ( 


P2-P, 


.^VV 


Thus,  Equations  (87),  (90)  and  (9I)  represent  4 ordinary  differential 

■ -j 

equations  in  4 unknowns,  6,,  6_,  9 and  4 . However,  due  to  the  present 

1/00 

assumption  of  a single  component  perfect  gas,  62  ” “<5|  necessarily  follows; 
accordingly.  Equation  (87)  need  be  solved  only  for  6|,  reducing  the  system 
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to  3 ode's,  a standard  4th  order  Runge-Kutta  integration  method  has  been 
used  to  integrate  the  system,  given  appropriate  initial  conditions. 

6.4  Numerical  Solution  - The  conditions  corresponding  to  the 
muzzle  blast  data  of  Schmidt  and  Shear  shown  in  Figure  (32)  were  selected 
to  test  the  model.  The  initial  conditions  were  determined  from  a laminar 
solution  of  the  heat  condition  problem: 


T = i ^''’2‘^1^  [ 1 + erf  y/(2/^)]  (93) 

Using  T = 0.995  (^2*^1^  '^1  gives  a laminar  solution: 

6,  = 3.65Arr'  (94) 

The  laminar  viscosity,  specific  heat  and  Prandtl  number  were  then  used  to 

estimate  a"  initial  value  of  i-j  . Consequently,  it  was  only  necessary  to 

select  arbitrarily  the  initial  ° value  of  0^  \ fortunately,  as  will  be 

shown,  the  growth  rate  of  6^  is  insensitive  to  arbitrarily  small  initial 

o 

values  of  0^  • 

The  transition  time  t at  which  the  turbulent  solution  should  be  initiated 
was  estimated  at  about  400  usee,  using  the  departure  of  the  experimental 
data  from  the  inviscid  contact  surface  solution,  shown  in  Figure  (32).  As 
demonstrated  in  Figure  (33)  the  resulting  solution  for  6j  is  in  reasonably 
good  agreement  with  the  experimental  observations.  (The  last  two  data  points 
represent  the  limit  of  the  field  of  view;  the  arrows  indicate  that  the  tur- 
bulent front  was  somewhere  beyond  this  distance.)  However,  by  increasing  the 
transition  time  to  500  Jsec.  the  solution  can  be  brought  into  virtually  per- 
fect agreement  with  the  data,  while  by  decreasing  it  to  about  200  usee,  the 
growth  of  6j  is  overestimated;  note  however  that  the  comparison  is  made  on 
a log-log  plot.  Fortunately,  the  dependence  on  assumed  transition  time  dimin- 
ishes with  increasing  time  and  vanishes  asymptotically. 


"71- 


TR  2i*6 


2 

It  can  be  seen  from  Figure  (3^)  that  the  RMS  temperature  fluctuation 
predicted  by  the  present  two-equation  TTE  model  asymptotes  to  about  10% 
of  the  mean  temperature  (T^+T2)/2.  Furthermore,  the  eddy  conductivity  K^, 
as  shown  in  Figure  (35),  appears  to  asymptotically  approach  a slowly  varying 
function  of  the  quantity  p(6^-62)  ^ suggests  that  a classical 

eddy-viscosity  or  eddy-conductivity  first-order  closure  model  may  be  adequate 
(asymptotically)  for  this  problem. 

This  numerical  example  and  comparison  with  experimental  data  is  regarded  as 
confirming  the  validity  of  the  present  model.  The  extension  to  a chemically 
reacting  mixture  could  proceed  along  the  following  line.  Invoking  a 
"turbulent  Lewis  number"  to  define  the  eddy  diffusivity  of  the  species,  it 
is  only  necessary  to  add  the  species  continuity  equation,  including  the 
diffusion  terms,  to  the  present  system.  A species  distribution  function, 
similar  to  Equations  (85)  or  (88)  must  be  defined,  e.g.. 


a.  “ o.  (1 -cos2irn) /2  + [a.  n + (1-n)a.  ] (l  + cos2irn)/2 
' 'o  'l  '2 


(95) 


The  species  continuity  equation  need  then  be  evaluated  only  at  y=0,  to 
determine  requiring  one  additional  ODE  to  be  solved  for  each  species 

in  the  reaction  model.  However,  even  with  the  relatively  simple  mechanism 
given  by  Equations  (6^a)  through  (6l4h) , some  eleven  species  and  accordingly 
eleven  additional  equations  are  involved.  Although  the  present  integral 
method  is  computationally  very  efficient  compared  to,  say,  a finite  differ- 
ence scheme,  it  was  not  considered  warranted  in  view  of  the  exploratory 
nature  of  the  study  and  the  absence  of  comparative  experimental  data,  to 
carry  out  any  nonequilibrium  chemistry  examples. 


J 
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VI  I . CONCLUSIONS 

Formation  of  an  unsteady  plume  and  blast  field  following  ignition  of  a rocket 
motor  has  been  described.  !t  is  pointed  out  that  the  supersonic  plume  which 
is  believed  to  form  along  the  ground  surface  initially  expands  (i.e.,  the 
Mach  disc  moves  outward)  as  the  blast  over-pressure  decays,  but  then  later 
collapses  (i.e.,  the  Mach  disc  falls  inward  toward  the  origin)  as  the  rate  of 
energy  deposition  at  the  ground  surface  drops  off.  The  rate  and  duration  of 
energy  deposition  have  been  shown  to  be  highly  dependent  on  turbulent  mixing 
in  the  steady  exhaust  stream  of  the  rocket  motor  (prior  to  impingement  on  the 
ground)  for  low  acceleraticn  (a  < 2g)  launches.  On  the  other  hand,  very  simple 
inviscid  flow  approximations  for  the  rate  and  duration  of  energy  deposition  ap- 
pear appropriate  for  very  high  acceleration  (a  > 20g)  launches.  The  net  effect 
is  to  reduce  the  dependence  of  energy  deposition  on  vehicle  acceleration  in 
the  range  of  intermediate  rates  (2g  < a < 20g) . This  quantitative  assessment 
of  viscous  effects  is,  of  course,  somewhat  dependent  on  the  particular  example 
considered  which  involved  a single  equivalent  nozzle  approximation  for  the  ex- 
haust stream  of  a multi-nozzle  engine  in  the  2x1(J^N  thrust  category. 

The  asymptotic  behavior  of  the  blast  field  (for  both  M5>>1  and  Mg-*!)  has 
been  examined,  and  the  present  numerical  solution  of  the  inviscid,  unsteady 
flow  between  the  blast  wave  and  Mach  disc  has  been  shown  to  be  in  excellent 
agreement  with  asymptotic  theory  in  both  limits.  The  latter  limit  (Mj-^l ) 
provides  the  more  stringent  test  of  the  numerics.  Correlations  of  numerical 
solutions  for  the  shock  layer  over-pressure  and  density  distributions  in 
terms  of  the  Vihitham  variables  are  presented,  and  may  in  point  of  fact,  be 
utilized  in  lieu  of  the  numerical  solution  for  Ms<1.1, 

The  formulation  of  the  numerical  method  of  solution  for  the  unsteady  layers 
between  the  blast  wave  and  contact  surface  and  between  the  contact  surface 
and  Mach  disc  has  been  presented.  Stability  of  the  method  subject  to  wave 
propagation  and  chemical  rate  limitations  has  been  discussed.  Interestingly, 
when  both  limitations  are  about  equal,  the  maximum  permissible  time  step  is 
about  one-half  of  that  due  to  either  individual  limitation;  only  when  the 
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two  are  of  distinctly  different  magnitudes  does  one  or  the  other  dominate 
the  stability  of  the  system.  An  example  is  presented  to  characterize  the 
pressure  and  temperature  field.  The  persistence  of  a major  temperature  Jump 
across  the  contact  surface  is  emphasized,  as  well  as  the  absence  of  a velo- 
city jump.  The  possibility  that  the  temperature  jump,  or  more  specifically 
the  associated  density  jump,  could  produce  a Taylor  instability  of  the  con- 
tact surface  is  pointed  out. 

Extension  of  the  fluid  dynamic  model  to  include  nonequilibrium  chemical 
reactions  has  been  discussed.  Several  numerical  examples  are  presented  to 
illustrate  the  thermochemical  characteristics  of  the  flow  field.  The 
assumption  of  a freeze  out  of  chemical  reactions  from  the  nozzle  exit  to 
the  ground  surface  plume  Mach  disc  has  been  shown  to  be  suspect.  In  the 
absence  of  mixing  across  the  contact  surface,  the  reactions  may  be  expected 
to  equilibrate  on  the  time  scale  of  the  subject  problem.  However,  vigorous 
afterburning  of  the  propellant  gases  may  be  anticipated  as  turbulent  mixing 
entrains  shock-heated  air  into  the  propellant  gas. cloud.  The  afterburning 
(and  attendant  optical  flash)  will,  therefore,  probably  be  mixing-controlled, 
but  the  possibility  of  chemical  non-equi librium  in  the  mixing  layer  in  its 
early  stages  should  be  considered. 

An  unsteady  turbulent  mixing  layer  analysis  has  been  developed  to  describe 
conduction  and  diffusion  across  the  contact  surface.  The  temperature  grad- 
ient is  assumed  to  be  the  driving  mechanism  for  the  production  of  turbulent 
thermal  energy  (RMS  temperature  fluctuations).  A two-equation  closure  model 
for  the  turbulent  thermal  energy  (TTE)  is  proposed,  by  analogy  with  Spaulding's 
two-equation  model  for  the  turbulent  kinetic  energy  (TKE)  equation.  A two 
strip  integral  method  solution  for  the  mixing  layer  has  been  developed. 

Growth  of  the  mixing  layer  thickness  is  found  to  be  in  excellent  agreement 
with  experimental  observations  of  a muzzle  blast  flow  field.  Extension  of 
the  method  to  a chemically  reacting  system  has  been  outlined. 

A complete  description  of  the  unsteady  plume  and  blast  field  evolution  follow- 
ing Ignition  and  (vertical)  1 junch  of  a rocket  has  been  presented.  The  meth- 
odology for  quantitative  analysis  of  the  blast  over-pressure,  afterburning 
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of  the  propellant  gases,  and  associated  optical  flash,  subject  to  the  princi- 
pal assumption  of  spherical  symmetry,  has  been  developed.  However,  the 
considered  processes  represent  only  the  initial  formation  of  the  propellant 
gas  cloud  at  the  ground  surface  following  lift-off.  The  ultimate  motion  and 
dispersion  of  this  propellant  gas  cloud  will  result  from  the  combined  effects 
of  buoyancy  and  winds,  which  have  not  been  considered  in  this  study. 
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